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ABSTRACT 


This  note  considers  the  problem  of  estimation  of  an  unknown  probability- 
density  function,  p(x),  and  the  derivative  of  the  logarithm  of  that  density 
function  , 


f(x) 


given  N  samples  from  an  ensemble  whose  probability  density  is  p(x).  Our 
principle  objective  is  to  estimate  f(x)  since  it  has  been  shown  that  the 
optimal  receiver  for  known  threshold  signals  in  additive  white  (but  possibly 
non- Gaus sian)  noise  consists  of  a  filter  matched  to  the  signal  preceded  by  a 
no-memory  device  whose  transfer  characteristic  is 


f(x)  | 


z 


x=  y 


where  y  is  the  input  to  the  device  and  z  the  output. 

Although  our  approach  is  primarily  motivated  toward  obtaining  a  good 
estimate  of  f(x),  the  method  and  results  would  appear  also  to  be  applicable 
to  finding  "good"  estimates  of  p(x).  The  "goodness  "  of  the  estimation  pro¬ 
cedure  is  investigated  theoretically  and  experimentally. 
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I. 


INTRODUCTION 


This  note  considers  the  problem  of  estimation  of  an  unknown 
probability  density  function,  p(x),  and  the  derivative  of  the  logarithm  of  the 
density  function 

f(x)  =  ^  [  In  p(x)  ]  (1.1) 

given  n  samples  from  an  ensemble  whose  probability  density  function  is  p(x). 
Our  principle  objective  is  to  find  a  "good"  estimate  f(x)  because  of  its  use 
in  communications  theory  applications.  For  example,  it  has  been  shown 
that  the  optimal  receiver  for  known  threshold  signals  in  white  (but  possibly 
non-Gaussian)  noise  is  a  filter  matched  to  the  signal  preceded  by  a  no¬ 
memory  device  whose  transfer  characteristic  is 

Z  =  f(x)  I  (1.2) 

x=y 

where  y  is  the  input  to  the  device  and  z  the  output.  Although  our  approach 
is  primarily  motivated  toward  obtaining  a  good  estimate  of  f(x),  the  method 
and  results  would  appear  to  be  applicable  to  obtaining  ’’good"  estimates  of 
p(x).  ' 

The  topic  of  obtaining  "good"  estimates  of  the  density  function  has  some 
relevance  for  our  particular  problem  and  so  we  first  review  some  approaches 
proposed  by  other  investigators.  Next,  we  present  our  estimation  procedure 
which  is  based  on  the  assumption  that  p(x)  is  smooth  enough  over  a  finite 
interval  ("the  smoothing  interval")  such  that  lnp(x)  can  be  accurately 
represented  by  a  power  series  with  a  finite  number  of  terms.  This  allows  us 
to  characterize  the  problem  of  finding  a  "good"  estimate  of  p(x)  and  f(x)  as 
one  of  estimating  the  power  series  coefficients. 

We  find  that  a  power  series  expansion  with  two  terms  gives  a  readily 
implemented  estimation  procedure  (based  on  a  maximum  likelihood  estimation) 
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An  approach  somewhat  similar  to  ours  has  been  suggested  by  Liden  for 
estimating  p(x)  in  connection  with  pattern  recognition  studies. 
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for  estimating  p(x)  and  f(x)  as  well  as  lower  bounds  of  the  variance  of  the 
estimates.  The  "goodness"  of  the  estimates  in  terms  of  bias  and  variance  is 
examined  experimentally  using  sample  functions  drawn  from  known  ensembles. 

Finally,  we  indicate  possible  use  of  this  procedure  in  communications 
theory  and  pattern  recognition  studies. 

II.  REVIEW  OF  PROBABILITY  DENSITY  FUNCTION  ESTIMATION 


Despite  the  importance  of  problems  related  to  density  function 
estimation,  there  is  comparatively  little  literature  regarding  the  estimation 
procedure  for  the  case  in  which  one  has  only  a  very  vague  idea  of  what  the 
form  of  the  density  function  might  be.  At  the  outset,  it  must  be  noted  that 
all  smoothing  operations  of  interest  assume  some  regularity  of  the  function 
being  estimated,  e.g.  the  existence  of  derivatives,  or  the  vanishing  of  high 
order  Fourier  coefficients  etc.  We  somewhat  arbitrarily  classify  prior 
work  we  are  aware  of  into  three  classes: 
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1.  "Distribution  independent 11  linear  smoothing  [Bartlett,  Parzen, 

g 

Rosenblatt  ]  :  These  papers  consider  the  class  of  estimates  of  p(x)  of  the  form 

+u(n)/2 

p  (x)  =  f  ,  v CO  (x-s)dP  (s)  (2.1) 

where  P^(x)  is  the  sample  distribution  function  for  n  samples.  If  co^(s)  is 
taken  to  be  suitably  well  behaved  even  function  of  s  such  that 

u(n)/2 

,  \  /o  |co  (s)  ds  =  1 
-u(n)/2  J  n 


it  can  be  shown  that  by  having  the  width,  u(n),  of  the  smoothing  interval  go 
to  zero  as  n-*®  such  that  nu(n)-*CD,  then  £  (x)  converges  to  p(x)  in  quadratic 


n 


mean  for  all  points  of  continuity  of  p(x)  with  the  variance  of  £  (x)  going  to 


,  ,i 

1  -1 

zero  as 

nu(n)  J 

L 

n 


*  7  s 

Parzen  shows  that  it  is  sufficient  to  require  that  go  (s)  =  K — - — r  where  K(y) 
.  j-i.-  n  u(n) 

meets  the  conditions  :  _ 

00 

1.  suP_co<v<oo  lK(y)  I <0°  2.1im  m  fyK(y)  J  =  0  3.  J  |K(y)  |dy<® 
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2.  "Partial  distribution  dependent"  linear  smoothing  [Whittle^]  : 

The  class  of  estimates  of  p(x)  considered  is  again  of  the  form 

u(n)/2 

p  (x)  =  f  t  ,/  go  (x-s)  d  P  (s).  (2.2) 

*nv  J-u(n)/2  nv  n  v  ' 

However,  in  this  case,  the  go  (s)  are  chosen  to  minimize  the  variance  of 

n 

p^(x)  over  a  given  family  of  density  functions  [an  example  of  this  would  be  to 

obtain  the  appropriate  Gon(s)  when  the  apriori  p(x)  is  taken  to  be  Gaussian 

with  unit  variance  and  a  mean  uniformly  distributed  between  +10  and  -10]  . 

For  the  classes  of  p(x)  considered,  it  is  found  that  the  variance  of  p  (x) 

-  1  n 

does  not  go  to  zero  faster  than  [n] 
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3.  "Distribution  dependent"  smoothing  [Cramer,  Parzen,  and 
many  others]  .  In  this  case,  one  assumes  a  model  of  p(x)  that  is  completely 
specified  except  for  the  values  of  a  finite  set  of  parameters.  The  values  of 
these  parameters  are  then  estimated  by  methods  such  as  maximum  likelihood 
estimation  or  the  method  of  moments.  In  a  case  where  sufficient  evidence 
exists  to  suggest  that  a  model  containing  a  manageable  set  of  parameters 
will  fit  the  data,  this  method  probably  offers  the  best  procedure.  However, 
in  many  cases  of  interest,  the  number  of  parameters  required  and  the  form 
of  equations  that  result  in  obtaining  parameter  estimates  make  this  method 
difficult  to  apply  in  practice. 

The  first  two  approaches  are  analogous  to  determination  of  a  suitable 
linear  estimator  of  a  continuous  waveform  given  a  received  waveform  whereas 
the  third  approach  is  analogous  to  determining  the  best  (e.g.  ,  maximum 
likelihood  and/or  minimum  variance)  nonlinear  estimator  of  the  continuous 
waveform.  For  example,  in  the  case  where  the  apriori  density  function  is 
taken  to  be  Gaussian  distributed  with  zero  mean,  the  third  approach  (using 
maximum  likelihood  estimation)  gives: 


where 


3 


A 


(2.4) 


CD 


We  note  that  the  estimator  of  equations  2.  3  and  2. 4  is  not  a  linear  functional 
of  the  sample  distribution  function  P^jx). 

For  our  particular  problem,  the  first  two  approaches  discussed  above 
have  the  liability  that  a  good  estimate  of  p(x)  may  not  yield  particularly 
good  estimates  of  f(x)  =  ^  ^Idx^X^  *  Furthermore,  it  is  not  exactly  clear  how 
one  would  operate  on  a  "good*1  estimate  of  p(x)  to  give  a  good  estimate  of 
f(x).  The  third  approach  gives  a  straightforward  procedure  for  estimating 
f(x);  but  as  we  have  indicated,  it  is  often  quite  difficult  to  find  a  manageable 
form  that  fits  the  data  over  the  required  range  of  x. 

HI.  ESTIMATION  PROCEDURE 


In  this  section,  we  present  the  assumptions  made  as  to  the 
regularity  of  p(x)  and  develop  in  detail  an  estimation  procedure. 


The  regularity  condition  that  we  impose  is  to  consider  the  logarithm  of 
the  density  function,  ln[p(x)]  ,  to  be  a  reasonably  smooth  function  of  x,  so 
that  in  the  vicinity  of  some  value  of  x,  say  x  =  x^,  the  function  ln[p(x)J  can 
be  expanded  in  a  power  series  in  x  -  x^  with  a  finite  number  of  power  series 
coefficients.  In  this  way,  estimating  the  value  of  p(x)  (and  parameters  such 
as  ~^[lnp(x)]  )  at  xq  reduces  to  estimating  the  power  series  coefficients  for 
the  expansion  around  x^ . 

A  model  is  assumed  for  the  distribution  of  the  sample  density  function, 


After  this  work  was  completed 

M 

proposed  the  model  lnp(x)  =  L 

i=  1 


I  became  aware  of  the  work  of  Liden^  who 
®  i^i(  x)  where  cp^(x)  represent  fixed  functions 


and  the  0.  are  estimated  from  the  data.  Liden  points  out  that  his  model  is 
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motivated  by  a  result  of  Dynkin  that  if  there  is  a  sufficient  statistic  of 
finite  dimension  M  that  characterizes  the  density  function  p(x),  then  p(x) 
must  be  of  the  form  above.  However,  Eiden  was  not  able  to  solve  the  resulting 
estimation  equations  for  any  cases  that  would  be  appropriate  to  our  problem. 
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p^(x),  conditioned  on  the  "true11  value  of  this  density  function  being  p(x). 

This  allows  us  to  obtain  the  likelihood  ratio  for  the  observed  sample  density- 
function  distribution  conditional  on  the  power  series  coefficients  having 
certain  values.  It  is  then  easy  to  obtain  maximum  likelihood  estimates  of 
the  power  series  coefficients  as  well  as  Cramer-Rao  bounds  on  the  variance 
of  the  estimates. 

From  the  assumption  that  ln[p(x)  ]  can  be  expanded  in  a  power  series 

of  finite  order  M  around  x  we  have: 

o 


lnp(x)  =  lnp(xj  +  d*1jxP^X>^  lx_  x  (X-XQ)  + 

o 


.  .  .  .  + 


dM  1  [  In  p(x )  1 


dx 


M-  1 


.  (x-v 

x  —  x  (M-l)  ! 
o 


M-  1 


(3.  1) 


M-  1 

l 

i  =  0 


a.  (x  )  (x-x  ) 
iv  o  o 


where  we  have  written  a^(xQ)  to  emphasize  that  the  power  series  coefficients 
will,  in  general,  be  a  function  of  xq.  Equivalently,  we  can  rewrite  (3.  1)  as: 

M-l 

p(x|a)  =  exp  [  ^  a.  (xq)  (x-x^1]  (3.2) 

i  =  1 


where  we  have  written  p(x  |a  )  to  indicate  that  p(x)  is  a  function  of  the  power 

series  coefficients  fa.}. 

1 J 

Next  we  consider  an  appropriate  probabilistic  model  for  the  distribution 
of  the  sample  density  function  PR(X)  conditional  on  p(x).  In  the  cases  of 
interest  to  us,  the  sample  density  function  p^(x)  can  be  obtained  from  the 
relationship. 


Pn(x)  = 


P  (x+-^-)  -  P  (x 
n _ 2 _ n 

Ax 


(3.3) 
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where  the  Ax  is  chosen  small  enough  so  that 


p(x  | a  ) 


,  Ax  ,  -v  r..  Ax 
P(x  +  —  |  a  )  -  P(x  -  —  |a  ) 

Ax 


3.4) 


where  P(x  |a  )  is  the  distribution  function  corresponding  to  p(x  |a  ).  The 

quantity  y(x)  =  np  (x)Ax  represents  the  number  of  samples  of  x(t)  whose 

Ax  Ax 

values  lie  in  the  interval  (x  -  -j-  ,  x+  — — ). 

Let  us  define  the  Smoothing  interval  M  as  the  closed  interval 
u(n) 


x 

o 


x 

o  2  J 


along  the  x-axis  in  which  equation  (3.1)  holds.  Next  we  divide  the  smoothing 

interval  into  Ns  non- over  lapping  intervals  (i.e.  ,  "bins")  such  that  u(n)/Ns 

equals  the  quantity  Ax  of  equations  3  and  4.  If  successive  values  of  x(t)  are 

independent,  the  joint  distribution  of  the  Ns  values  of  y(x)  corresponding  to 

the  various  bins  in  the  smoothing  interval  is  then  the  multinominal  distribution 

4 

with  Ns  +  1  exclusive  events  and  n  trials  [Cramer  (1948)]  .  The  likelihood 
ratio  for  the  observed  sample  density  function  values  Pn(x  )  in  the  smoothing 
interval  is  then: 


where  y 


A 


m 


1  =  l 


n p  (x  )Ax. 
rn  m 


n  ! 

'Ns 

n 

[  A x p  (x.  |a  )] 

yj" 

Lj=  i 

y.! 

J 

Ns 

Ns 

_  n 

-z 

1  -  z 

1=  1 

Ax  p  (x^  |a  ) 

JL  = 

(3.  5) 


The  maximum  likelihood  estimates  (ML  estimates  )of  the  a  are  obtained 
from  the  set  of  simultaneous  equations 


ST  {pr  (?„<•>  I-  !  }  =  °  to 


r  i  =  0,  1 _ M-l 


(3.6) 
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Substituting  (3.5)  into  (3.6),  we  obtain  the  set  of  equations 


Ns 

Z 

j=  i 


Ns 

A*  n  '  £=  X 


p(x|a*)  Ns  _ 

J  1  -  y  Ax  p(x ^  |a  )J 

1=  1 


Ax 


dp(x.  la 

_ _ JL_ 

da. 

l 


=  0 


(3.7; 


for  i  =  0,  1,  2,  .  .  .  M-  1. 

Equation  (3.7)  cannot  be  conveniently  solved  for  cases  of  interest.  However, 
for  the  usual  case  of  n  large  and 

Ns 

V  Axp(x^|a  )  <  <  1 

j 1=  1 


we  can  set 


Ns 


Ns 

1  -  Yj  AxP  (xjfc  la 

1=  1 


=  n 


(3.8) 


so  that  equation  (7)  becomes 


Ns 


Pn(xi) 

n  J  ■  -  1 


p(x-  la 
J  J 


Ax 


dp(x.  |a  ) 

da. 


=  0 


(3.  9) 


For  later  use  in  bounding  the  variance  of  our  estimates  of  the  a  ,  we 
note  that  (under  the  assumption  of  equation  8)  the  elements  of  F isher 's 


Information  Matrix  [Van  Trees  '( 1968)]  are 
*  E[ 


d£np[pn(.  |a  )]  d£np[pn(.  |a  )] 


da 


daT 


-  V  (n 


Ax) 


d p(x.  |a  )  dp(x.  |a 

-  M  -  J  - 

da  * 


(3.  10) 


da. 


i  p(Xj  |a  ) 


See  Appendix  A  for  a  discussion  of  this  assumption. 
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From  equations  (3.2)  and  (3.8),  we  find  that  the  ML  estimates  of  the  power 
series  coefficients  at  x=x^,  are  found  as  the  solution  to  the  set  of  equations 


z< 


x. 

J 


V 


i  r  =  ? 


M-  1 


a  (x. 
r  J 


V 


=  Bxj-xo)1pn(xj) 


(3.  11) 


where  i  =  0,1,2,....,  M-  1 .  If  M  is  allowed  to  be  greater  than  2  ,  this  set  of 
equations  cannot  be  solved  conveniently.  Thus  we  limit  the  investigation  to 
the  case  M=  2.  This  restraint  on  M  then  sets  a  limit  on  the  {  x j  |  to  be 
considered  in  obtaining  the  estimate  of  the  |a1|  at  x  =  x^:  we  consider  only 
those  |  xj }  which  lie  sufficiently  close  to  such  that  the  density  function 
can  be  adequately  described  by  the  power  series  expansion  of  (3.2)  with  M  =  2  . 


In  a  later  section,  we  indicate  the  results  of  some  experiments  as  to 
what  a  reasonable  maximum  range  of  {  x j  |  around  xq,  i.  e.  ,  size  of  u(n),  is 
for  some  well-known  density  functions.  In  cases  where  n  is  very  large,  one 
might  well  choose  a  u(n)  that  depends  on  n. 


In  Appendix  B,  we  show  that  if  the  density  function  p(x)  and  its  first 
two  derivatives  are  continuous  and  if  we  shrink  the  smoothing  interval  width, 
u(n),  down  as  n  1  where  0  <  a  ^  <  1/3,  then  the  estimates  of  p(x)  and 
f(x)  =  d/dx  np(x)J  using  M  =  2  will  be  unbiased  with  variances  that 
as  n-* 00 . 

For  the  case  M  =  2,  it  is  useful  to  define  p^  =  ea°,  so  that  p(x  |a  ) 
becomes 


p(x.  [p  a,)  =  p  ea^XJ  X°^  (3.12) 

r  j  ,ro »  1  ro 

and  p  ,  the  ML  estimate  of  pQ,  can  be  interpreted  as  the  ML  estimate  of 
p(x)  at  x  =  xq.  Substituting  (3.  12)  into  (3.  11)  and  carrying  out  the  algebra, 
we  obtain 
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(3.  13) 


(3.  14) 


where  the  function  g(.  )  is  defined  by  the  relationship 

a  (x. -x  ) 


g(a)  =  ^ 


E,xj-xo|e  ' J  ° 


j 


a  (x.-x  ) 
J  o 


(3.  15) 


and  is  the  ML  estimate  of  f(x)  at  x=x^.  The  function  of  g(n)  of  equation  (3.  15) 

must  be  obtained  numerically,  but  is  a  function  only  of  the  jx.-x^J  and  not 

x  alone.  Thus,  if  for  all  x  ,  we  keep  the  set  ix.-x  \  constant,  e(n)  and 
o  o  <  j  o  » 

its  inverse  need  to  be  computed  only  once  as  we  compute  a^  and  a^  for  all 

v 

In  Fig.  1  we  show  g(n)  for  u(n)  =  1*0,  i.e.,  |x.|  covering  the  x-axis 
interval  (x^  -  .5,  x  +  .5).  This  particular  u(n)  has  no  particular  significance 
until  the  relation  of  x  to  the  statistics  of  x(t)  is  known.  One  x-axis  scale 
for  which  Fig.  1  is  particularly  relevant  for  has  X  =  value  of  x(t)  in  root- 
mean  square  units  from  the  mean;  i.e.,  X  =  [x(t)  -  x(t)]  /&  ^t)' 


IV.  BOUNDS  ON  THE  ESTIMATION  ERRORS 

Next  we  consider  the  Cramer-Rao  bounds  on  the  variance  of  p 

r  o 

and  a^  assuming  that  p^  and  a^  are  unbiased.  Using  equations  (3.2)  and 

^Appendix  B  shows  that  by  having  u(n)  and  Ax-*0  in  an  appropriate  fashion  as 
n-*00,  that  p  and  a^  are  asymptotically  unbiased.  For  cases  where  n  is  large, 
but  finite,  it  is  difficult  to  analytically  establish  whether  or  not  pQ  and  a^  are 
unbiased  and/or  the  degree  to  which  the  Cramer-Rao  bounds  are  satisfied  with 
equality  in  equations  (4.  1)  and  (4.3).  Thus,  we  argue  that  the  utility  of  these 
bounds  can  best  be  established  by  a  consideration  of  how  well  they  predict  the 
observed  deviations  for  the  sample  functions  drawn  from  known  density  functions. 
In  a  later  section  we  present  some  simulation  results  that  indicate  that  the 
bounds  of  equations  (4.  1)  and  (4.  3)  with  equality  are  fairly  good. 
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4  9 

(3.  10),  and  the  well  known  '  expression  for  the  Cramer-Rao  bound  with  two 
estimated  parameters,  we  obtain: 


Var  (p0>  2  J  3 — ‘.‘j  j 

00  11  J01  10 


(4.  la) 
(4.  lb) 


(nAx)  ^eai^xj"Xo)J 


1- 


V/  \  ai(x.-x  ) 

)  (x.-x  )e  1  j  o 
Zj'  jo  j 


~t  2 


V  ai(x.-x  )  .2  a,(x.-x  ) 

2j  e  1  J  o  >J(xj"Xo)  e  1  J  ° 


LJ 


|  nDPo  eai(Xj"X°)Ax  |  j  1 


(4.  lc) 


E,  >2  ai(x.-x  )  v  ai(x.-x 

(Xj-xo)  e  j  o'  2,e  1V  j  o 


where  the  were  defined  by  equation  (3.  10).  It  should  be  noted  that  the 

first  term  in  the  denominator  of  equation  (4.  lc)  is  the  expected  number  of 
samples  of  x(t)  whose  values  lie  in  the  "smoothing  interval",  i.e.  ,  the  x 
axis  interval  spanned  by  the  Because  of  the  normalization  implicit  in 

the  ML  estimate  of  p^  (eq.  3.  14),  this  first  term  is  numerically  equal  to  the 
number  of  samples  of  x(t)  whose  values  lie  in  the  "smoothing  interval".  The 
second  factor 


Vo<al>  =  1  - 


E.  .  a,(x.-x  ) 

(x.-xo)e  1  j  o 


[2<x.-*o>2eal<*j-V]  |)>al'Xj-Xo>] 


(4.2) 


in  (4.  lc)  can  be  interpreted  as  a  multiplicative  factor  by  which  the  number  of 

samples  is  reduced  because  of  the  slope  (i.e.  ,  non -zero  a  ^ )  of  in  (  p(x)] 

near  x  =  x  . 

o 
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Similarly,  the  Cramer-Rao  bound  on  the  variance  of  a^  is  given  by 


,,  -  ,  .  _ 00 

Var  [ a  j ]  _ 

J11J00J01  10 


(4.3; 


The  first  bracketed  term  in  Eq.  4.  3  is  numerically  equal  to  the  average  number  of 
samples  of  x(t)  in  the  "smoothing  interval",  while  the  second  term 

V  /  a  l  (x.  -x  ) 

6<Xj-xo)  e  J  o' 

J _ 

Sai(x.-x  ) 
e  AV  J  °; 

j 

can  be  interpreted  as  a  multiplicative  factor  by  which  one  multiplies  the 
number  of  samples  in  the  smoothing  interval. 

In  Figs.  2  and  3,  we  plot  the  multiplicative  factors  depending  on  a^, 

V^(a^)  and  V^(a^),  that  arose  in  the  expressions  for  the  variance  of  p^  and 
a^  respectively  for  a  symetrical  "smoothing  interval".  Again  we  point  out 
that  the  Cramer-Rao  bounds  represent  lower  bounds  so  that  the  usefulness 
of  equations  4.  1  -  4.4  is,  at  this  point,  inconclusive.  In  a  later  section, 
we  will  present  results  that  indicate  that  the  bounds  established  here  are 
meaningful. 

V.  APPLICATION  TO  DETERMINISTIC  SAMPLE  DENSITY  FUNCTION 

As  indicated  earlier,  the  key  assumption  that  j£n[f(x)]  can  be  adequately 
represented  by  a  power  series  of  low  order  was  primarily  motivated  by  the 

The  smoothing  interval  width  of  0.  5  used  for  figures  2  and  3  might  be 
appropriate  in  the  case  where  x  represents  rms  units  from  the  mean. 


vi<ai 


)  * 


S,  x  a,(x.-x  ) 

(x.-xo)e  1  j  O 


V  a[(x.-x  ) 

Le  j  ° 


ii 


observation  that  the  logarithm  of  many  common  density  functions  is  indeed  a 
smooth  function  and  by  the  fact  that  this  assumption  leads  to  computationally 
feasible  results.  Thus  we  argue  that  an  appropriate  criteria  for  the  utility 
of  this  approach  is  how  well  it  does  at  estimating  pQ  and  a^  for  various 
density  functions. 

In  Figs.  4-7  we  show  the  results  of  applying  the  estimation  procedure 
to  a  (deterministic)  sample  density  function  whose  form  is 

Pjxi)  =  ~~~~  exP  [--5x;2]  (5.1) 

J  yz  n  J 

for  various  sizes  of  the  "smoothing  interval".  In  Figs.  4  and  6,  the  "x" 
represents  p(x)  at  various  x  values  while  the  smooth  line  represents  the 
m.l.e.  pQ(x).  In  Figs.  5  and  7,  the  light  line  represents  the  actual  slope  of 
£n[p(x)] 


_d_ 

dx 


[  jfcnp(x)]  lx=x 

o 


-  X 

o 


(5.2) 


while  the  dark  line  represents  a^(x),  the  ML  estimate  of  the  slope  of  |np(x). 
By  examining  the  differences  between  the  "true"  values  of  the  estimated 
parameter  and  the  ML  estimates  of  these  parameters  as  a  function  of 
"smoothing  interval"  size  for  various  density  functions,  one  can  obtain  a 
measure  of  the  bias  to  be  expected  using  the  procedure  outlined  here. 

The  results  for  the  Gaussian  fit  is  also  of  interest  because  e  can 

be  expanded  about  any  point  x=xq: 

2/,  x  ^  (x-x  )^ 

in  e  =  -  xq(x-xo) - ^ -  •  (5.3) 

We  note  that  had  we  included  the  (x-xq)^  term  in  our  power  series  expansion 

for  inp(x),  we  could  estimate  the  p(x)  |  _  and  -jr-  Anp(x)|  _  exactly  for 

x—  x  dx  x  ~  x 

o  o 
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a  sample  density  function  of  the  form  of  Eq.  (5.  1).  The  deviations  of  the 
M  L  estimates  from  the  "true  "  values  of  Eqs .  (5.1)  and  (5.2)  as  shown 
in  Figs.  4-7  thus  give  an  indication  of  the  error  one  might  expect  to  arise 
from  using  a  power  series  expansion  of  lnp(x)  consisting  of  only  two  terms. 
From  equation  (5.  3),  we  see  that  the  error  gets  larger  as  |X“XQ  I  increases; 
the  particular  smoothing  intervals  used  for  Figs.  4-7  were  chosen  larger 
than  one  might  use  in  practice  [particularly  near  the  mode  of  p(x)]  . 

VI.  APPLICATION  TO  NOISE  SAMPLE  FUNCTIONS 

Although  in  Appendix  B  we  establish  that  the  estimates  of  p(x) 
and  f(x)  are  asymptotically  unbiased  and  normally  distributed  with  variances 
given  by  equations  (4.  1)  -  (4.4),  it  is  extremely  difficult  to  analytically 
establish  any  useful  results  in  the  small  or  even  medium  sample  cases.  Thus, 
we  now  present  some  results  of  applying  the  estimation  procedure  to  sample 
density  functions  drawn  from  a  variety  of  noise  processes.  In  each  case  we 
show: 

1.  a  plot  of  the  sample  density  function  Pn(x)  and  the  estimate  PQ(X) 
of  the  density  function 

2.  a  plot  of  the  estimate  PQ(X)  and  the  "true  ",  i.e.  ,  ensemble  density 
function  p(x)  together  with  the  Cramer -Rao  variance  bounds  on  pQ(x).  The 
variance  bound  curves  plotted  correspond  to  p(x)  ±  3  */var  [  p^(x)]  where 
var(p^)  is  given  by  equation  (4.  1). 

3.  a  plot  of  the  estimate ,  a^(x),  of  f(x)  and  the  "true",  i.e., 
ensemble  f(x)  together  with  the  Cramer-Rao  variance  bounds  on  a^(x).  The 
variance  bound  curves  plotted  correspond  to  f(x)  ±  3  ,/var  [a  j(x)]  where 
var(a^)  is  given  by  equation  (4.  3). 

4.  a  plot  of  the  expected  number  of  samples  in  the  smoothing  interval 
centered  at  x  as  a  function  of  x. 

5.  a  plot  of  the  sample  exceedance  probability,  1  -  FJx),  and  the 
ensemble  exceedance  probability,  1  -  F(x),  to  give  some  information  as  to  how 
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"typical  11  the  particular  sample  function  used  was. 

Figures  8-12  show  the  results  for  n=  50,000  and  u(n)  =  1.0  rms  units 
from  the  mean  for  a  random  number  simulation  of  a  noise  process  with  a 
probability  density  function  of  the  form  p(x)  =  A  for  x  ^  0.  This 

particular  density  function  is  of  interest  because  the  density  function  in  any 
region  of  x  >  0,  will  be  fit  exactly  by  our  two  term  power  series  representa¬ 
tion  of  ln[p(x)]  .  We  note  that  both  the  estimates  are  significantly  different 
from  the  true  values  in  the  region  of  width  1. 0  corresponding  to  x  =  0.  This 
arises  because  the  step  in  p(x)  at  x  =  0  makes  the  regularity  assumption 
invalid.  Over  the  range  of  x  for  which  our  regularity  condition  does  hold, 
the  estimates  and  "true  "  values  agree  to  within  the  variance  bounds. 

Figures  13-17  show  the  results  for  n  =  50,000  and  u(n)  =  1.0  rms  units 
from  the  mean  for  a  random  number  simulation  of  a  Gaussian  noise  process. 

In  this  case,  there  are  no  steps  or  other  abrupt  changes  in  the  ensemble 
density  function  and  so  our  regularity  assumption  is  approximately  valid  over 
the  entire  range  of  x.  In  Figs.  13-17,  we  see  that  the  estimates  agree  with 
the  "true"  values  to  within  the  variance  bounds. 

Figures  18-22  show  the  results  for  n  =  50,000  and  u(n)  =  1.0  rms  units 
from  the  mean  for  a  random  number  simulation  of  a  Pareto  process  with 
density  function  p(x)  =  3 /x  for  x  >  1.  As  in  the  case  of  the  exponential 
process  of  Figs.  8-12  we  see  the  effect  of  the  step  in  p(x)  at  x  =  1  is  to 

a  a 

cause  the  estimates  PQ(X)  and  a^(x)  to  differ  significantly  from  the  "true" 
values  over  an  interval  of  width  u(n)  centered  at  the  discontinuity.  There  is 
also  some  indication  in  Figs.  18  and  19  that  the  estimators  are  slightly  biased 
(in  the  sense  of  giving  estimates  that  consistently  are  slightly  outside  the 
variance  bounds  in  one  direction).  However,  this  bias  (if  it  indeed  exists) 
would  appear  to  be  quite  small. 

VII.  DISCUSSION 

In  this  section,  we  discuss  a  minor  extension  of  our  estimation 
procedure  as  well  as  some  possible  uses  of  the  procedure  in  communications 
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and  pattern  recognition  problems.  In  dealing  with  density  functions  of  low 
frequency  electromagnetic  noise  (which  is  known  to  be  non-Gaussian),  we 
have  found  it  helpful  to  use  a  small  smoothing  interval  to  estimate  data  near 
the  mode  (since  the  density  function  changes  very  rapidly  near  the  mode) 
and  a  much  larger  smoothing  interval  for  estimation  on  the  tail  of  the  density 
function  (since  there  is  a  far  smaller  number  of  samples  per  unit  of  smoothing 
interval  length  out  on  the  tails).  Since  the  estimation  procedure  is  done 
separately  for  each  bin  of  the  density  function,  it  is  easy  to  make  the  smooth¬ 
ing  interval  width  a  function  of  distance  from  the  mode  [although  function 
g  ^  (a)  of  equation  (3.  15)  must  be  recomputed  every  time  the  smoothing 
interval  width  is  changed]  . 

We  have  already  mentioned  the  utility  of  f(x)  in  receiver  design.  It 
can  be  shown  (as  will  be  done  in  a  subsequent  report)  that  the  function  f(x) 
is  also  of  use  in  estimation  in  non-Gaussian  noise  as  well  as  in  obtaining 
estimates  of  the  error  performance  improvement  that  can  be  obtained  using 
an  optimal  receiver  in  non-Gaussian  noise.  Finally,  we  note  that  our 
estimation  procedure  might  also  be  of  use  in  pattern  recognition  studies  of 
the  type  considered  by  Liden.  ^ 
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APPENDIX  A 


In  this  appendix,  we  discuss  briefly  the  rationale  behind  the  assumption 

that 


Ns 

n  -  i  y* 
_ 


Ns 

1  -  Yj  AxP(X£  la  ) 
JL=  1 


-  n 


(Al) 


used  in  obtaining  a  solvable  set  of  maximum  likelihood  equations.  Our 

4 

discussion  makes  heavy  use  of  three  results  in  Cramer  namely  that 


e[Y£]  =  nAxP(xJa)  ( A2 ) 

E  {[y£  "  E(y  l  )1  2}  =  n(Axp(x^|a)]  [  1  -  ax  p(x^  |a  )]  (A3) 

and 

E{[y£-E(y£)]  [ym-E<ym)]}=  -n(Ax)2  p(x£  |a)p(xm|a) 

for  H  rh  m. 

We  note  that  equation  (Al)  can  be  written  in  the  form 


"“-SV '  2(yi'npil 

>-2p. 


where  p  -  ixp(x,  la). 

I  l' 

Equation  (A4)  clearly  reduces  to 


(A4) 
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(A5) 


n[1  -l  E(lHnp*.>1 

”  ‘-Ep, 


Let  us  examine  the  statistics  of 


R  ^ 


i  t(yf-npt' 

n  1  •  Ep« 


(A6) 


From  (A Z),  it  is  clear  that 


E[R]  =  0. 

Using  (A3)  and  (A4),  we  have 


Var [ R]  =  — 


1 


£ 


{EEI'V,-"P£>2i} 


Ns  Ns 


E  ,  E  ,  Ei<yi-"pjt»  <ym-npm>i 

o =  1  m=  1  * 


l 

1 


n  (1  -ZP,)  Li 

l 


— 

Ns 

Ns 

EnP£(1-P£> 

i 

-z 

Z  np£p 

m=  1  1  ' 

l 

Ns 

X=  1 

cp£<*-p£»-; 

5P*(- 

pjj+  z ,  ■ 

x  m=  1 

nd  -  >,  ?/ 


[■5pi(£ip“,  +  2pJ 


Thus  we  have 


18 


Ns 

<Z  .V 

Var  [  R]  =  — -  (A7 ) 

n(l-Jp  ) 
i 

_ (Probability  of  x(t)  in  smoothing  interval) 

(number  of  samples )(probability  of  x(t)  not  in  the 

smoothing  interval) 

Typically,  the  probability  of  x(t)  in  the  smoothing  interval  is  less  than 
0.2,  so  that  the  standard  deviation  of  R  is 

aRs—  .  (A8) 

For  cases  where  the  procedure  outlined  here  might  be  used,  we  would 
expect  n  to  be  at  least  100,  for  which 


aR  £  -05  .  (A9) 

Thus,  setting 

n(l  +  R)  =  n  (A10) 

should  be  a  reasonable  approximation  for  a  "worst  case"  of  interest,  and 
represents  an  approximation  that  improves  as  we  increase  n. 
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APPENDIX  B 


Asymptotic  Behavior  of 


/N 


and 


1 


In  this  appendix,  we  show  that  if: 

1.  the  ensemble  density  function  p(x)  and  its  first  two  derivatives 
are  continuous  and  finite 

2.  the  width,  u(n),  of  the  smoothing  interval  goes  to  zero  propor- 
tional  to  n  1  as  n  (the  number  of  independent  samples  of  x) 
goes  to  oo  (where  0  < 


We  can  choose  the  bin  width,  £x,  so  that  the  estimates  of  p(x)  and 
f(x)  =  ^  [ln  p(x)]  will  be  asymptotically  unbiased  with 


-  1 


•  *  1  -  fY  i 

1.  the  variance  of  the  estimate  of  p(x)  going  to  zero  as  (n  1 ) 

2.  the  variance  of  the  estimate  of  [lnp(x)]  going  to  zero  as 

,  1  •  3of  i  v  -  1 
(n  A)  . 


From  the  conditions  on  p(x),  we  can  represent  P(x)  in  the  smoothing 
interval  by  the  power  series 

(x“Xo)2 

p(x)  =  P(xQ)  +  a1p(xo)(x-xQ)  +  p"(xo) - - -  (Bl) 


where 


&  P'(X°)  d  „  ,  . 

ai  -  rtTT  =  di  f ln p(x))  !x=x 
r  o  o 


(B2) 


The  proof  goes  as  follows:  we  will  first  show  that  a  random  variable 

closely  related  to  the  righthand  side  of  equation  3.  13  is  asymptotically 

normal  with  a  mean  equal  to  the  true  value  and  a  variance  that  goes  to  zero 
1  3  1 

as  (n  1)  .  From  this,  we  will  show  that  the  estimate  &  ^  is  therefore 

asymptotically  normal  with  a  mean  value  equal  to  and  a  variance  that  goes 
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I  1 

to  zero  as  (n  al)  .  Finally,  we  will  show  that  the  estimate  is 
asymptotically  normal  with  a  mean  equal  to  p^  =  p(x^)  and  a  variance  that 
goes  to  zero  as  (n^ 

First  we  consider  the  statistic 

Ns 

Y  (x.-x  )y . 

U  .  i  o,y j 
1=1  J _ _ 

Ns 

2  yj 

j=  i 


where : 


R,  = 


12 

2,  \ 
Lu  (n). 


(B3) 


Ns  =  number  of  bins  in  u(n)  = 
A 


x 


X. 

J 


=  number  of  values  of  x(t)  that  lie  in  the  bin  corresponding 

to  x. 

J 

=  Xo  +  (1  '  T]  AX 


S  S  - Y  (x.-x  )y. 

„u3(n)  f.l  3  °  1 

a  ,  Ns 
S  £  —Lf  V  y. 
y  n  u(n)  U  1  7j 


(B4) 


(B5) 


Ax  =  bin  width  =  k,  n  a2  (0  <  <  1)  .  (B6) 

Let  us  define 


pj = 

Ax  p(x^  |a  ) 

j  =  1,2, .. 

.  ,  Ns 

(B7) 

qj = 

1  -Pj 

II 

LI 

Ns 

l-f  P.  = 

f=  i  J 

xQ+u(n)/2 

1  -  J* 

1 XQ  -u(n) /2 

p(x  |a  )dx 

(B8) 

^L  = 

1  -  PL 

(B9) 

21 


It  is  clear  that  as  n 


E(y-j)  =  npj  =  nAxp(x.)-.nAxp(xQ)  =  k^n1  "Q’2  p(xQ)  .  (BIO) 

-♦  5 

Thus,  the  y  meet  the  conditions  of  Gnedenko's  "local  limit  theorem" 

from  which  we  obtain 


p(y|a)  = 


■  1  2 
Ns  -  —  q.  z . 

n  ^ _ li_ 

1  Jtn  P: 


J 


1  2 
2  qLZL 


(Bll) 


where  z.  =  (y.-np.)/y^q~  . 


A  NS 

=  (n  -  V  y .  - 

r=  i  J 


L  .  ’)  nPL'/V“KIj'iL 

J=  1  J 

As  n-»co,  we  have  from  (B6)  through  (BIO): 


)/J™ 


Pt  qi 


q-  1 


qL  =  “*  0  • 


Thus 


N 

p(y  |a  )  =  n  e 

j=  1 


.1  (yi-np/ 


np. 


(B 12) 


,/2  77  np. 


5, 

from  which  we  conclude  the  j  y .  |  are  asymptotically  normal  and  independent 

An  interesting  side  note  that  indicates  the  rate  at  which  the  y.  become 
independent  is  the  behavior  of  correlation  coefficient  between  y.  and  y^ 

(j  4  k)  as  n->oo 

E{[yj-E(y-)]  [yk-E(yk)]j  np^pk _ 

0  jk  rV^  Nr  \  vTl  Nr  M  172  (n  r,  n  r^  n  1  ^ 


[Var  (y.)  Var  (yk)]  [n“p.q^pkqk] 


- -i 


_  -  Axp(x  )  =  k  n  “2  p(x  ) 
qjqk  H  0  2  o' 
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with 


E(y.)  =  np. 


(B  1  3) 


Var  (y  )  =  np  .  (B14) 

Thus,  the  sums  S  and  S  are  asymptotically  normal  with  means  and 
variances  that  are  to  be  calculated. 


Let  us  now  compute  the  mean  and  variance  of  and  respectively. 

From  the  body  of  the  report  and  Appendix  A,  we  recall  that  the  y  have  a 
multinomial  distribution  corresponding  to  n  repetitions  and  Ns  +  1  exclusive 
events.  We  will  work  out  E(S^^)  in  detail  to  indicate  the  approach  and  then 
state  the  results  for  Var(S  ),  E(S  )  and  Var(S  ).  Let  us  define  x.  =  x. -x  . 


Then: 


xy 


J  J  o 


^  E< v  =  Z  s5E(yj)  3  s  *J  npj 


Ns 

Z 

j=i 


Ns 

Z 

j=  i 


Ns 


=  Yi  *j  AxP(Xj  |a  )] 


i=  1 


(B 15) 


Ns  p"(x  )  2 

=  V  x.nix  [p(x  )  +  a ,  p(x  )(x.-x  )+ - —  (x.-x  )  1 

U  j  1  rv  o  lr  o'  jo  2  j  o  1 

j=l 


which,  as  n-.ro,  becomes 


u(nl 


P"(xo>  2, 


ilih(ll]-E(Sxy)  =  nf  u(n)  *  [P(x0)  +  a1p(xo)(x)+-T-^(x)  ]  dx 


linl 


=  n  [-®ip(x  )  +&£  a  p(x  )  +-^  p"(x  )\U*  ] 

2  0  3  1  0  4  o  u(n) 


=>  nalP(xo)  -L^- 


(B16) 
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so  that 


E(S  )  =  a.  p(x  )  . 

xy  1  rv  o 

Similarly,  in  the  limit  as  u(n)-+  0,  it  can  be  shown  that 


(B  17) 


E(Sy)  =  p  (xo) 

P(Xo} 

Var<Sy>  =  ^) 


Var(S  )  = 


p(xQ) 


xy  5 .  x 

J  nu  (n) 

Substituting  u(n)  =  k^n  °*1  into  (B19)  and  (B20),  we  find  that  as 


(B 18) 
(B  19) 

(B20) 


n_«°° 


Var  (S  ) 

y 


(B21) 


12  p(x  ) 

Var (S  )  -  -  7  °. 

xy  ki3nl-3or1 


(B22) 


From  (B21)  and  (B22),  we  conclude  that  as  n-»  oo ,  S  converges  in 

y  4 

probability  to  the  constant  p(xQ).  Thus,  from  Cramer  (pages  254-55)  ,  we 
conclude  that  the  density  function  of  (S  /S  )  *s  asymptotically  normal  with 


E(S  ) 
xy 


E^ac)  =  ,  . 

Sy  P(Xo} 


=  a . 


Var  (-ft)  = 


Var(S^)  [12/k1ip(xo)] 


1  -3q;  . 
n  1 


(B23) 


(B24) 


y  [p(x0)J 

Next  we  wish  to  consider  the  statistical  behavior  of  the  estimate  of  a 
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1  f  u2(n)  Sxyl 
L  12  s  J 


1 


(B25) 
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where 


As  n-*cof 


q(r>)  = 


Ns 

Z<Y 

_  1=  i  j 


x  )e  J 


q(r»)  =  i 
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Ns 
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r> 
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J  Ax 
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J  ° 

Ns 

v* 

T7^j-x0) 

L  e 

Ax 

j=  i 

u(n) 

2 
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)  J*  x 

e^X  dx 

2 

(B26) 


»x 


u(n) 
r  2  e71"'  dx 

J_uin) 

2 


(B27) 


Let  be  the  largest  (finite)  value  of  ^  for  which  we  wish  to  know  the 
value  of  q(f)).  Then  as  n-*» 


l»?*l  5  =  1^0  I  (_2')n’ftl  ° 

so  that  we  can  make  the  substitution 


(B28) 


=  1  + 


V 


X 


and  thus  obtain 


(B29) 


so  that  for 


/  v  u  (n) 
q(r?)  =  12  *1 


2.  ,  S  2/  > 

u_ (n)  xy  u  (n) 

12  S  *  12  ^o 


(B30) 
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we  have 


•,  -  *»  ¥  *vo  ■ 

y  y 


(B32) 


Now  from  our  previous  results  for  the  ratio  S  /S  ,  we  realize  that  by 


xy'  y 


choosing  77  >  a  ^ ,  we  can  ass  ert  that  for  0  <  a  j  <"3* 
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Pr 


^O^  = 


0 


(B33) 


n-»  00 


Thus,  is  asymptotically  normal  with  mean  a^  and  variance 

2  _  _ 12 

o-  L  3  l-3rv 

1  k,  p(x  )  n  1 

1  r  v  o 


(B34) 
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[expected  number  of  samples  of  x(t)  in  the  smoothing] 
interval  u(n) 

Finally,  we  consider  the  statistical  behavior  of  our  estimate  of  p(x^) 
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(B35) 


(B36) 


From  our  results  regarding  the  asymptotic  behavior  of  a^,  we  know  that  for 
any  Aj  >  a 
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CD 


Pr(a"1  >  Aj) 


0  as  n-» 


Thus,  from  (B5),  we  have 
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expected  number  of  samples  of  x(t) 
in  the  smoothing  interval  u(n) 
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Fig.  1.  G(cv)  versus  cv. 
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Fig.  2.  Density  function  variance  multiplicative 

factor  V  (a,)  versus  a,. 

o'  1  1 
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Fig.  3.  [in  p(x)]  variance  multiplicative  factor  V^(a^)  versus  a^. 
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PROBABILITY  DENSITY  FUNCTION 


Fig.  4.  Results  of  density  function  estimation  for 
deterministic  sample  density  function. 
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Fig.  5.  Results  of  [in  p(x)]  estimation  for  deterministic 
sample  density  function. 
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PROBABILITY  DENSITY  FUNCTION 


Fig.  6.  Results  of  density  function  estimation  for 
deterministic  sample  density  function. 
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Fig.  7.  Results  of  [in  p(x)]  estimation  for  deterministic 
sample  density  function. 


34 


PROBABILITY  DENSITY  FUNCTION 


RMS  UNITS  FROM  THE  MEAN 


Fig.  8.  Sample  estimated  density  functions  for 
exponential  ensemble  density  function. 
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RMS  UNITS  FROM  THE  MEAN 


Fig.  9.  Comparison  of  estimated  and  ensemble  density  functions 
for  exponential  ensemble  density  function. 
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RMS  UNITS  FROM  THE  MEAN 


Fig.  10.  Comparison  of  estimated  and  ensemble  ^  [In  p(x)]  for 
exponential  ensemble  density  function. 
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EXPECTED  NUMBER  OF  EVENTS  IN  SMOOTHING  INTERVAL 


RMS  UNITS  FROM  THE  MEAN 


Fig.  11.  Expected  number  of  samples  in  smoothing  interval 
for  exponential  ensemble  density  function. 
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EXCEEDANCE  PROBABILITY 


RMS  UNITS  FROM  THE  MEAN 


Fig.  12.  Comparison  of  sample  and  ensemble  exceedance 
probabilities  for  exponential  density  function. 
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PROBABILITY  DENSITY  FUNCTION 


RMS  UNITS  FROM  THE  MEAN 


Fig.  13.  Sample  and  estimated  density  functions  for 
Gaussian  ensemble  density  function. 
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PROBABILITY  DENSITY  FUNCTION 
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RMS  UNITS  FROM  THE  MEAN 


Fig.  14.  Comparison  of  estimated  and  ensemble  density 
functions  for  Gaussian  ensemble  density  function. 
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RMS  UNITS  FROM  THE  MEAN 


Fig.  15.  Comparison  of  sample  and  ensemble  [in  p(x)]  for 
Gaussian  ensemble  density  function. 
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EXPECTED  NUMBER  OF  EVENTS  IN  SMOOTHING  INTERVAL 


RMS  UNITS  FROM  THE  MEAN 


Fig,  16.  Expected  number  of  samples  in  smoothing  interval 
for  Gaussian  ensemble  density  function. 
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EXCEEDANCE  PROBABILITY 


RMS  UNITS  FROM  THE  MEAN 


Fig.  17,  Comparison  of  sample  and  ensemble  exceedance 
probabilities  for  Gaussian  ensemble  density  function. 
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PROBABILITY  DENSITY  FUNCTION 


RMS  UNITS  FROM  THE  MEAN 


Fig.  18.  Sample  and  estimated  density  functions  for 
Pareto  ensemble  density  function. 
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PROBABILITY  DENSITY  FUNCTION 


RMS  UNITS  FROM  THE  MEAN 


Fig.  19.  Comparison  of  estimated  and  ensemble  density 
functions  for  Pareto  ensemble  density  function. 
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Fig.  20.  Comparison  of  estimated  and  ensemble  [in  p(x)]  for 
Pareto  ensemble  density  function. 
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EXPECTED  NUMBER  OF  EVENTS  IN  SMOOTHING  INTERVAL 


RMS  UNITS  FROM  THE  MEAN 


Fig.  21.  Expected  number  of  samples  in  smoothing  interval 
for  Pareto  ensemble  density  function. 
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EXCEEDANCE  PROBABILITY 


Fig.  22.  Comparison  of  sample  and  ensemble  exceedance 
probabilities  for  Pareto  ensemble  density  function. 
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